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Simulation of a turbulent flame in a channel 

By G. Bruneaux 1 , K. Akselvoll 2 , T. Poinsot 3 AND J. H. Ferziger 2 


The interaction between turbulent premixed flames and channel walls is studied. 
Combustion is represented by a simple irreversible reaction with a large activation 
temperature. Feedback to the flowfield is suppressed by invoking a constant density 
assumption. The effect of wall distance on local and global flame structure is inves- 
tigated. Quenching distances and maximum wall heat fluxes computed in laminar 
cases are compared to DNS results. It is found that quenching distances decrease 
and maximum heat fluxes increase relative to laminar flame values. It is shown 
that these effects are due to large coherent structures which push flame elements 
towards the wall. The effect of wadi strain is studied in flame-wall interaction in 
a stagnation line flow; this is used to explain the DNS results. It is also shown 
that ‘remarkable’ flame events are produced by interaction with a horseshoe vortex: 
burnt gases are pushed towards the wall at high speed and induce quenching and 
high wall heat fluxes while fresh gases are expelled from the wall region and form 
finger-like structures. Effects of the wall on flame surface density are investigated, 
and a simple model for flame-wall interaction is proposed; its predictions compare 
well with the DNS results. 


1* Introduction 

The interaction of a turbulent premixed flame with a wall is quite complex. First 
the flame is strongly influenced by the presence of the wall; which limits flame wrin- 
kling and may cause the flame front to quench. Moreover, the flame has a significant 
effect on the flow in the vicinity of the wall: viscosity is greatly increased in the 
burnt gases, inhibiting turbulence. At the same time, flame elements approaching 
the wall increase the heat flux to as much as 1 MW/m 2 in practical situations. For 
these reasons, modeling flame-wall interactions in turbulent flows is an important 
issue (Amsden et al 1985, Clendening et al. 1981, Lu et al 1990). Models which 
try to predict these phenomena are available (Jennings 1992, Poinsot et al. 1993). 
However, little fundamental information is available so model building is a difficult 
exercise. An additional problem is that experiments are difficult to perform because 
the interesting phenomena occur very close to walls (typically less than 1 mm). 

Our objective is to explore the flame-wall interaction mechanisms using three- 
dimensional direct numerical simulation (DNS). Two-dimensional variable- density 
simulations were performed in 1992 (Poinsot and Haworth, 1992) and led to a model 
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used in piston engines (Poinsot et al. 1993). One of the main difficulties was the 
lack of a statistically stationary turbulent flow field. In the present study, a constant 
density turbulent channel flow was used. This has advantages and drawbacks: (1) 
the turbulence characteristics in the channel flow are well known and stationary 
which allows easy computation and (2), a constant density approximation has to 
be used, prohibiting feedback of the flame effects to the flow. However, this is a 
cost-effective approach. 


2. Numerical method and configuration 

In this study we extended the three dimensional DNS channel flow code written by 
Akselvoll & Moin (1993) to take reaction into account. Temperature and fuel mass 
fraction are treated as passive scalars and do not affect the flow. The flow solver 
has not been modified and is independent of the solver for the chemical species. 

2.1 Basic equations 

The flow solver solves the Navier- Stokes equations for an incompressible, constant 
viscosity flow: 


dpiii _ dp dp UjUi t d 2 ui 

~W = ~dTi dTj + **' dsf 


dxii 

dxi 


= 0 


( 2 ) 


The reaction solver solves the energy and species conservation equations, which 
allow convection, diffusion, and reaction effects: 


dpc p T dpc p uiT _ d 

dt dxi dxi 

dpY F dpujYp 
dt dxi 



(3) 

(4) 


The superscript (~) refers to physical variables; absence of a superscript indicates 
a dimensionless variable. _ _ 

We assume that pD == pD\D* where D* = ( T/T\) b and A = AiA* where A* = D *. 
The subscript 1 refers to the fresh gases, and the subscript 2 refers to the burnt 
gases. 

The reaction is represented by a simple one-step mechanism, corresponding, for 
example, to lean combustion in which fuel is the limiting factor in determining the 
reaction rate (Williams, 1985). The reaction rate is expressed as: 


wr = p YfB exp 



( 5 ) 


where B is the pre-exponential factor and T a is the activation temperature. 
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The equations are nondimensionalized using the following dimensional quantities: 
u r , the friction velocity at the wall, A, the channel half width, 7\, the temperature 
in the fresh gases, T 2 , the temperature in the hot gases, Y\, the fuel mass fraction 
in the fresh gases. 

Physical and dimensionless variables are related in the following way: 

U = U T U , t = ht/u T , p = fm* T p, X = hx, T = (f - 7\)/(t 2 - f,J, Y F = Y X Y F 
giving the set of dimensionless equations: 


dui _ dp dujui 1 d 2 u y 

dt dxi dxj ^ Re dx 2 



( 6 ) 

(7) 


dr du^T _ _i d_ / ar\ . 

dt + dxi RePr 5x, \ + WR 

OYf dufr = 1 d ( D *dY F \ _ . 

dxi RePrLe dx y \ dx y ) WR 


( 8 ) 

(9) 


where the Reynolds number is Re = the Prandtl number is Pr = and 

the Lewis number is Le = — . The reaction rate expression can be reduced to 

pepD 1 

(Williams, 1985): 

w R = DdYpexp (- ^^(1-T) ) (10) 

where a is the temperature factor a = (T 2 — Ti)/72, ft is the reduced activation en- 
ergy j3 = aT a /T 2 , and Da is the reduced pre-exponential factor Da = S^exp ^ . 


2.2 Numerical implementation 

These equations are solved in a Cartesian coordinate system using a second-order 
finite difference scheme. All terms are treated explicitly except the diffusive terms in 
the wall-normal direction in the momentum equation, which are treated implicitly. 
The time discretization is second-order Adams-Bashforth for the explicit terms and 
second order Crank-Nicolson for the implicit terms. The pressure is used to correct 
the velocity field so that it satisfies the continuity equation; this requires a Poisson 
solver. 

The flow and flame structure are computed on different meshes. The mesh 
spacing for the flow needs to be small enough near the wall to resolve the vis- 
cous sublayer, typically A = 0.1, but in the center of the channel, A = 8 
is sufficient (Kim et al. 1987); here the superscript ( + ) denotes wall units. Be- 
cause the structure of the turbulence is elongated in the streamwise direction, 
Ax + = 35 and A = 5. The mesh distribution along the y axis is given by 

A tanh(o^(l— 0— l)/(Ny-l»)^ Yr , . , , r 

yj = (^1 * — L J where a is a stretching parameter. Large 
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values of a distribute more points near the wall. NY is the number of points and 
Yi is the size of the box in the y direction. 

The flame computation is best done on a uniform mesh in all directions in order 
to resolve the flame equally well everywhere in the computational domain. 

For this reason, different meshes are used in the computation, one for the velocity 
and pressure, and another for the temperature and fuel mass fraction. The velocities 
are interpolated from the flow mesh to the reaction mesh. Three-dimensional linear 
interpolation based on data at the eight corners of the smallest box (of the flow 
mesh) surrounding the reaction grid point is used. This interpolation procedure 
was tested by running laminar flame-wall interaction cases in two dimensions and 
does not introduce additional error because the velocity grid near the wall is much 
denser than the temperature grid. 

2.3 Initial and boundary conditions 

The walls are no-slip and isothermal. The flow, temperature and mass fraction 
fields are periodic in the x and z directions. 

The initial conditions for the flow are obtained by running the flow solver until 
stabilized (in the statistical sense) values of the velocities and pressure are obtained. 
The temperature and mass fraction are introduced at t = 0 as two back-to-back 
one-dimensional laminar flames propagating towards the walls. 

3. Code validation and computation of reference flows 

S.l Computation of non- reacting turbulent channel flow 

A first calculation was made without a flame to validate the flow solver. The 
flow field is initialized from a random field, and is run until the values of velocities 
and pressure stabilize. Mean quantities are calculated by averaging in the x and z 
directions and time. 

The configuration is the minimal channel flow with Re = 180, and uses a stretch- 
ing factor a = 2.9. The dimensions of the domain are Xi = 3.14159, Yi = 2., and 
Zl = 0,908, with NX = 18, NY = 130 and NZ = 34. The results are consistent 
with well-known results of channel calculations (Kim et al 1987). Fig. 1 shows 
the mean velocity profile along with the log law, and Fig. 2 shows the profiles of 
turbulent velocity components compared to the results of Kim et aL 1987. The dis- 
crepancies are due to the fact that we have performed a minimal channel simulation 
while Kim et ai did a full channel simulation. 

S.2 Quenching of laminar flames on walls in stagnant flow 

One-dimensional calculations were performed to validate the reaction solver. The 
temperature profile of a constant viscosity flame is calculated using an analytical- 
numerical approach (Rutland 1989). Then this profile is used as an initial condition 
to compute a flame having variable transport properties with the reaction code, 
in which the flame is stabilized by prescribing a uniform inlet velocity equal to 
the flame speed. The temperature profile obtained is then used to initialize the 
flame-wall interaction calculation (laminar or turbulent). 
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FIGURE 1. Profiles of mean velocity compared to the laminar and log laws. 
linear law, log law, mean velocity. 

We first performed a one-dimensional calculation corresponding to head-on quench- 
ing in stagnant fluid. The flame propagates normal to the wall, fresh gases axe 
trapped between the flame and the wall, and the fluid velocity is zero everywhere. 
The flame consumes reactant as it moves towards the wall. When the flame-wall 
distance 6 reaches its minimum, the wall heat flux $ is maximum; the consumption 
rate s c decreases to zero exponentially thereafter, as shown in Fig. 3. The phenom- 
ena occurring in this interaction have been discussed previously (Adamczyk and 
Lavoie 1978, Carrier et al. 1979, Wichman and Bruneaux 1994). The flame-wall 
distance is non-dimensionalized by a typical flame thickness d = A /(pCpS ® ) to form a 
Peclet number Pe = 6 /d and the wall heat flux is non-dimensionalized by the flame 
power P = pc p s*f (T 2 — Ti) to produce a reduced heat flux <f> = $/P. At quenching 
the minimum Peclet number is Pe = 3.68 and the maximum reduced wall heat flux 
is <j> = 0.56. These values are different from the previous results of Poinsot et al. 
(1993) because the assumptions are different; in particular, they allowed variable 
density and had a different Prandtl number (Pr = 0.75). 
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FIGURE 2. Profiles of fluctuating velocities. DNS minimum channel, 

— Kim et al 1987. 

Table I. Fixed parameters for DNS of turbulent channel. 


Re 

Le 

Pr 

b 

a ft 

Da 

s°iu T 

6?/h 


180. 

i. 

0.5 

1 . 

0.75 6. 

80.4 

0.36 

0.15 

5.0 


We investigated the influence of grid resolution on the laminar flame speed (mass 
consumption rate) for a stabilized flame and on the maximum wall heat flux during 
flame-wall interaction. The results are shown in Fig. 4. For the maximum wall 
heat flux we also compared first and second order treatment of the wall boundary 
condition ( ] =0. The results show that 65 points in the half-channel and 

\ d * / wall 

a second order scheme at the wall suffice. 

3.8 Quenching of laminar flames on walls in stagnation line flow 
We also performed one-dimensional calculation of flame interacting with a wall 
in a stagnation line flow. A flow field corresponding to a stagnation line flow with 
= — r*, T, the strain rate is non-dimensionalized by the inverse characteristic 
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time/flame time 

FIGURE 3. Time evolution of wall heat flux, flame speed (consumption rate), and 
flame-wall distance for a laminar flame-wall interaction in a stagnant flow, 

flame time scale to produce a reduced strain rate 7 = r/(sj /£{*). A typical 

result for strained flame-wall interaction is presented in Fig, 5, for 7 = 5 . At 
first, the flame adjusts to the flow and the flame speed decreases to a stable value. 
Then the flame begins to interact with the wall. Because the flame is convected 
towards the wall, the interaction is faster than in the stagnant flow and produces a 
higher maximum wall heat flux and a smaller minimum flame wall-distance. In the 
stagnant case, the interaction lasts about 6 flame times, while for a strained flame 
with 7 = 1, it lasts 3.5 flame times, and for 7 = 5, 3 flame times. Fig. 6 shows the 
effect of the strain rate 7 on the maximum wall heat flux <f> and on the minimum 
flame-wall distance Pe. 


4. Results for turbulent flame wall interaction 

4-1 Evolution of global quantities during interaction 

Fig. 7 presents the evolution of the mean fuel mass fraction, the turbulent flame 
speed, and the total turbulent flame surface for a typical periodic run together with 
the laminar values. The parameters are those of the non-reacting flow and the 
laminar flame. 

Early in the simulation, the flames are not yet wrinkled and are located near 
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FIGURE 4. Influence of grid resolution and order of accuracy on laminar flame 
speed and maximum wall heat flux, x x first order, x x second order. 


the center of the channel where the turbulence is weak. After t/tf > 3 where 
tj = d[s{ = Ai /pCpis®) 2 , d is a typical flame thickness and s the laminar flame 
speed, wrinkling increases as does the consumption rate of reactants, producing 
a maximum turbulent flame speed 1.5 times the laminar flame speed. Then, at 
t/tf = 10, the flames begin to interact with the walls and the consumption rate is 
reduced due to lack of reactants. The turbulent flame speed is never constant. 

Fig. 8 presents the surface with T = .85 and the reaction rate in the turbulent 
flame at t/tf = 9. Quenching is observed, and a finger- like structure is also present. 
The mechanism of formation of this structure will be explained below. 

J.2 A correlation between local wall strain rate and flame quenching 

For practical applications the most important quantity is the maximum wall 
heat flux. Fig. 9 displays the variations of the minimum flame wall distance and 
the maximum wall heat flux (normalized by laminar quantities) with time. Fig. 9 
also displays the effect of grid resolution on these quantities: a second calculation 
was performed with twice the number of points in the y direction. No effect of 
grid resolution is seen. Clearly, the turbulent flame comes closer to the wall and 
produces higher heat fluxes than the unstrained laminar flame. This differs from 
results obtained previously and appears to be due to the structure of the turbulence. 
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FIGURE 5. Time evolution of wall heat flux, flame speed and flame- wall distance 

in a laminar flame- wall interaction in a stagnation line flow with 7 = 5. flame 

speed/laminar flame speed, — wall heat flux/flame power (</>), flame-wall 

distance (Pe). 

In the simulations of Poinsot et al (1993), the turbulence was two-dimensional 
so there was no small scale structure near the walls. In the present case, the 
typical quenching distance 8q is larger than the viscous sublayer thickness (typically 
8q ~ 28) and turbulent structures modify the structure of the flame near the wall. 

Another way of presenting this phenomenon is to look at trajectories in a (Peclet 
number-heat flux) diagram. Fig. 10 presents such trajectories for the laminar flame, 
for the turbulent flame at t/tf = 7.3, and for two strained laminar flames. The en- 
velop of the turbulent results lie close to the trajectory of the laminar strained flame 
with 7 = 5, indicating that the effect is primarily due to the strain. This was fur- 
ther checked by computing the strain rate statistics shown in Fig. 11. Because the 
strain rate is not constant in the turbulent case, we plotted the normal component 
of the velocity at the flame location divided by the flame-wall distance. We will see 
in the next section how these large strain rates are created. 

4.3 The importance of flow structures 

Near-wall coherent structures have a strong effect on the flame. In the case 
presented here, we find the interaction of a horseshoe vortex with the flame is quite 
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Laminar Peclet number in stagnant flow] 
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FIGURE 6 . Effects of wall normal strain 7 on flame- wall quenching of laminar 
flames. ■ wall heat flux/flame power (<f>), □ flame wall distance (Pe). 



FIGURE 7. Evolution of global quantities during turbulent flame wall interaction. 

(a) laminar, turbulent, (b) laminar flame speed, turbulent 

flame speed, -□ - total flame surface/initial flame surface. 
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region of low fuel 



Figure 8. Snapshot of isosurface of temperature T = 0.85, along with reaction 
rate isolines and fuel mass fraction field, t/tf =9. 

important. It produces the following events: 

(1) The horseshoe vortex pushes burnt gases towards the wall and leads to flame 
quenching with small Peclet number and large heat flux. 

(2) At the same time, the other side of the horseshoe vortex pushes fresh gas 
away from the wall, leading to the formation of an unburnt gas tongue; tongues 
similar to this one have been seen in experiments. 

Fig. 12 shows an instantaneous picture of a one-legged horseshoe vortex wrapping 
a flame and leading to quenching at the wall and ejection of fresh gas. 
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FIGURE 9. Time evolution of minimum flame-wall distance, minimum, mean, and 
maximum wall heat flux in the turbulent calculation. Comparison with a higher 
resolution calculation is also shown. NY — 130, — NY = 260. 

5. A model for the quenched interface density 

5.1 Mean quantities 

Since periodic boundary conditions were used, averaging may be performed in the 
two directions parallel to the wall (x and z ) at each instant. Quantities computed 
using conventional averaging include the mean fuel mass fraction Kf, the mean 
temperature T, and the mean reaction rate a?. It is convenient to replace the 
mean reaction rate a; by an equivalent reactive flame surface density £/* defined by 
£# = Zj/s*}. Profiles of these quantities are plotted in Fig. 13. 

We also estimated the density £ = (£') of interface between fresh and burnt 
gases. £' is the local surface to volume ratio, calculated where the surface (defined 
as the isosurface with reduced temperature 0.85) is approximated using the angle 
between the local temperature gradient and a coordinate direction jltutland 1989^ 
In the absence of quenching or strain, £ is related to u; by u; = £s{* or £# = £ 
(we used Lewis number unity for this flame. Near the wall, part of this interface 
is quenched so u? < £ 3 °. We will characterize quenching by the quenched fraction 

defined by Q = 

Early in the simulation ( t/tf = 1.8, Fig. 13), the flame is fax from the wall, no 
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flame- wall distance (Pe) 

FIGURE 10. Trajectory in the Pe-wall heat flux plane of the turbulent flame at 
t/tf = 7.3, along with trajectories for laminar flame- wall interaction in stagnant 
flow and stagnation line flow, o turbulent t/tf = 7.3, laminar. 

quenching takes place, and the profile of the interface density E matches the profile 
of the normalized reaction rate E/j. The noise in Fig. 13 is due to the fact that 
E is computed by estimating surface area while the mean reaction rate w and the 
reactive interface density Ej* are computed using conventional averages. The mean 
fuel mass fraction at the wall is still the initial value and the burnt gases occupy 
only a small fraction of the channel. The quenched fraction Q is zero everywhere; 
the wall heat flux is also essentially zero. 

Later (t/tf = 12.8, Fig. 13), the flame brush starts to interact with the wall 
and quenching takes place. This increases the interface density relative to the 
reactive interface density and thus the quenched fraction near the wall. The mean 
temperature gradient is not zero at the wall, indicating that the mean wall heat 
flux is no longer zero (see Fig. 9). In addition, the mean fuel mass fraction at the 
wall starts to decrease. 

Finally (t/tf = 14.7, Fig. 13), most of the fresh gases in the channel have been 
consumed and the interface density is much larger than the reactive surface density. 
Most of the interface is quenched. Very little fuel is available and the wall heat flux 
is large. 
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FIGURE 11. Trajectories in Pe- (equivalent) normal wall strain rate diagram at 
t/tf = 7.3. 


5,2 A model for the quenched interface fraction 

It is well known that interaction between flames and walls or, more generally, the 
behavior of non-adiabatic flames may be characterized in terms of enthalpy loss Lh 
(Williams 1985, Wichman and Bruneaux 1994) defined by 

L h = l-(»+T) (11) 

In an adiabatic premixed flame with unity Lewis number, Lh is zero everywhere. 
When the flame is non-adiabatic (as near walls), Lh increases, which indicates that 
quenching is possible. This is true for turbulent flames if we assume that heat and 
species diffuse at the same rate (turbulent Lewis number equal to unity). 

A simple model for the quenched interface fraction Q may be derived by assuming 
that Q is proportional to the enthalpy loss Lh times the interface density: 


Qmodti = Lh * 2 * Lq ( 12 ) 

where Lq is a multiple of the laminar flame thickness 6® , For the present, we assume 
Lq = l(wy. 
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FIGURE 12. Snapshot of isosurface of low pressure (which marks a horseshoe 
vortex) together with the velocity and fuel mass fraction fields at tjtf = 7.3. 

Despite its simplicity, this model retains much of the physics of flame quench- 
ing: quenching occurs only where flame surface is present and the flame has lost 
significant enthalpy. 

This model was tested against DNS results and the results are given in Fig. 14. 
The agreement between the modeled value Q mo del and the DNS value Q is good. 
The model predicts both the spatial extent of the quenching as well as its magnitude. 
Only late in the simulation ( t/tf = 14.7, Fig. 14) does the model underpredict the 
extent of quenching and then only in the region far from the wall. 

Conclusions 

Direct numerical simulations of flame-wall interactions have been performed using 
a three-dimensional channel flow code and a constant density reaction solver. Non- 
reacting turbulent flow and laminar reaction in stagnant and stagnation line flow 
were used to validate the code. In turbulent flows the wall heat fluxes are much 
higher than in the laminar stagnant flame-wall interaction. This is due to the 
turbulence which convects flame elements towards the wall, inducing high heat fluxes 
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FIGURE 13. Mean values of temperature, mean fuel mass fraction, reactive surface, 
and interface density, a): t/tf = 1.8; b): t/tf — 12.8; c): t/tf = 14.7. 
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Figure 14. Comparison between modeled and DNS-measured quenched fractions. 
-□-DNS mesured quenched fraction, -■-Model for quenched fraction, a): t/tf — 9; 
b): t/tf = 12.8; c) t/tf = 14.7. 
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to the wall. The turbulent flame was compared to a flame in a stagnation line flow, 
leading to the conclusion that high wall heat fluxes are due to high normal strain. 
In the turbulent case, the high normal strain is generated by horseshoe vortices 
which push flame elements towards the wall while fresh gases are convected away 
from the wall, forming finger-like structures. A model for the quenched fraction of 
interface was proposed and compares well to the DNS results, despite its simplicity. 
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